Phase Bubbles and Spatiotemporal Chaos in Granular Patterns 



Sung Joon MooniJ] M. D. Shattuck\0 C. Bizon^, Daniel I. Goldman\ J. B. Swift\ and Harry L. Swinney^EI 

^ Center for Nonlinear Dynamics and Department of Physics, University of Texas, Austin, TX 78712 
^ Colorado Research Associates Division, Northwest Research Associates, Inc. 3380 Mitchell Lane, Boulder, CO 80301 

(February 1, 2008) 

We use inelastic hard sphere molecular dynamics simulations and laboratory experiments to 
study patterns in vertically oscillated granular layers. The simulations and experiments reveal that 
phase bubbles spontaneously nucleate in the patterns when the container acceleration amplitude 
exceeds a critical value, about 7g, where the pattern is approximately hexagonal, oscillating at one- 
fourth the driving frequency (//4). A phase bubble is a localized region that oscillates with a phase 
opposite (differing by tt) to that of the surrounding pattern; a localized phase shift is often called an 
arching in studies of two-dimensional systems. The simulations show that the formation of phase 
bubbles is triggered by undulation at the bottom of the layer on a large length scale compared to 
the wavelength of the pattern. Once formed, a phase bubble shrinks as if it had a surface tension, 
and disappears in tens to hundreds of cycles. We find that there is an oscillatory momentum 
transfer across a kink, and this shrinking is caused by a net collisional momentum inward across 
the boundary enclosing the bubble. At increasing acceleration amplitudes, the patterns evolve into 
randomly moving labyrinthian kinks (spatiotemporal chaos). We observe in the simulations that 
//3 and //6 subharmonic patterns emerge as primary instabilities, but that they are unstable to the 
undulation of the layer. Our experiments confirm the existence of transient //3 and //6 patterns. 

PACS number(s): 45.70.Qj, 47.54. -hr 



I. INTRODUCTION 

Spatiotemporal chaos, where the physical variables 
vary in time and space in a seemingly random way, may 
arise when a spatially extended system is driven far from 
the primary instability [Q. In some systems, spatiotem- 
poral chaos is understood as a complicated evolution of 
the amplitude field, and the behavior can be described 
in terms of the dynamics of this field In some other 
systems, the dynamics of defects plays a primary role 
in the transition to spatiotemporal chaos; cases include 
Rayleigh-Benard convection electrohydrodynamic 
convection chemical patterns and Faraday insta- 
bilities . Several mathematical models have been pro- 
posed to describe the transition to spatiotemporal chaos 
in spatially extended physical systems, including ampli- 
tude chaos 0, phase turbulence j|] , defect-mediated tur- 
bulence pO|p| , and invasive defects but the under- 
standing is still far from complete. 

In this paper, we study patterns around the transi- 
tion from the //4 subharmonic hexagonal pattern to spa- 
tiotemporal chaos in granular layers with large aspect 
ratio {L/H ^ 10, where L is characteristic horizontal 
size of the layer, and H is the depth of the layer). The 
layers are subject to a sinusoidal oscillation in the direc- 
tion of gravity. The oscillation is characterized by two 
nondimcnsional control parameters, T = An^f'^A/g and 
/* = fy/H/g, where A is the amplitude of the oscilla- 
tion, / = 1 /T is the frequency of the oscillation, T is the 
period of the oscillation, and g is the acceleration due to 
gravity. We also define the nondimcnsional depth of the 
layer, N = H/a, where a is the diameter of the parti- 
cle. Various subharmonic standing wave patterns have 



been observed as a function of F and /* |13|; however, 
the transition to spatiotemporal chaos has not previously 
been investigated. In this study, we show that the transi- 
tion to spatiotemporal chaos in granular patterns is due 
to the intrinsic dynamics of oscillated granular layers; a 
large lengthscale undulation of the layer and an oscilla- 
tory momentum transfer across a kink. 

The rest of the paper is organized as follows. Section 
II presents the methods in the simulation and the ex- 
periment. Kinks, phase bubbles, and randomly moving 
labyrinths are described in Section III. In Section IV a 
large lengthscale undulation of the layer and its relation 
to the nucleation of a phase bubble is described. Section 
V discusses how a phase bubble shrinks and disappears. 
In Section VI, the prediction and observation of tran- 
sient //3 and //6 patterns are presented, and the paper 
is concluded in Section VII. 



II. METHODS 

A. Numerical simulation 

In the absence of well- validated macroscopic governing 
equations for vertically oscillated granular layers, current 
theoretical investigation proceeds at a more basic level, 
that of individual particles. Bizon et al. developed 
an event-driven inelastic hard sphere molecular dynamics 
simulation of this system, by implementing the collision 
operators in Ref. ||l^. This collision model conserves 
both linear and angular momentum, but allows energy 
to be dissipated through inelastic collisions and surface 
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friction. The normal coefRcient of restitution e{vn) de- 
pends on the magnitude of the normal component of 
relative colliding velocity w„ = (vi — V2) • fi2, where 
ri2 = (ri — r2)/|ri — ^2]- The accurate form of e{vn) is 
not known; it was assumed that e{vn) = 1 — Bvn^^ for 
Vn less than a crossover velocity Vc, and e(w„) = Co oth- 
erwise. The value of B was set to make e(w„) continuous 
at Vn = Vc- Here we use Vc — ^Jgo and = 0.7. The 
simulation results are not sensitive to the form of e(u„) 
for Vn < Vc- The tangential impulse is given by a coeffi- 
cient of friction \i times the normal impulse, with a cutoff 
corresponding to the crossover from a sliding contact to a 
rolling contact. The crossover ratio of the relative surface 
velocity after the collision to that of before the collision, 
—/3c, was set to —0.35 as suggested in Ref. and was 
set to 0.5. The values of these parameters were chosen 
to fit the wavelength of the pattern obtained in the ex- 
periments with lead particles, for three different control 
parameter sets. These fitting parameters reproduced the 
observed patterns quantitatively throughout the control 
parameter space [Q. The collisions between the grains 
and the container were treated in the same way as the 
collisions between grains. The mass of the container was 
assumed to be infinitely large compared to that of the 
granular layer. 

We performed three types of simulations: (1) 2D, (2) 
quasi-2D, and (3) 3D. We simulated a 2D or quasi-2D 
layer as vertical cross-section of a 3D layer (Figs. 3 and 
10). A quasi-2D layer is a 3D layer whose dimension 
in one direction is short enough lOtr) to be homoge- 
neous in that direction; these simulations run much faster 
than a fully 3D simulations yet yield the same statisti- 
cal information (Figs. 7, 8, and 11(a)). We performed 
simulations of 3D layers of square-shape with horizontal 
periodic boundary condition (Fig. 11(b)) and of 3D lay- 
ers of cylindrical shape with side wall, when we compare 
with the experiments (Figs. 4, 9). 



B. Experiment 

Experiments were conducted with vertically oscillated 
layers of granular material consisting of spherical bronze 
particles of mean diameter 165 [im (Spherical lead parti- 
cles of diameter 165 \xm were used only for Fig. 1). The 
nondimensional depth of the layer, iV, was in the range 
of 5 to 15, and the aspect ratio LjH ranged from 40 to 
150. Both circular and rectangular containers with var- 
ious sizes were used in the experiments. The container 
was mounted on an electro-magnetic shaker, and it oscil- 
lated sinusoidally in the direction of gravity with a single 
frequency, in the range 10 - 150 Hz. The value of F varied 
from to 14. The container was evacuated to a pressure 
of 4 Pa to reduce the role of interstitial gas. The con- 
tainer was encircled by a ring of LEDs, and the images 
were taken by a digital camera mounted above the con- 
tainer. A more detailed description of the experimental 



apparatus is found in Ref. [Q. 



III. PATTERNS AROUND A TRANSITION TO 
SPATIOTEMPORAL CHAOS 

The phase diagram of the patterns in oscillated gran- 
ular layers has been reported previously ||l^,|l^,|l^ . We 
present a new phase diagram in Fig. 1 that shows the 
details above the //4 subharmonic hexagonal pattern 
regime, which is the focus of the current study. 



A. Single inelastic ball model and temporal 
dynamics of the layer 

Much of the dynamics of the patterns in this sys- 
tem can be understood from the single inelastic ball 
model |l^Jl^. In this Section we review the results of 
this model in the //2 and //4 patterns regime. The 
single ball model is a one-dimensional model of the os- 
cillated granular layer, which approximates the center of 
mass of the layer as a completely inelastic ball (e = 0) 
on an oscillating plate. 

For F > 1, the magnitude of the acceleration of the 
container exceeds g during a fraction of the cycle, so that 
the layer loses contact with the container when the plate's 
acceleration becomes — and then the layer makes a free 
flight until colliding with the container later. In the //2 
square/stripe pattern regime, the flight time of the layer 
is a fraction of the oscillation period T, and the layer 
leaves and hits the container every cycle (Fig. 2(a)). In 
this regime, the magnitude of the acceleration of the con- 
tainer at the collision is less than g (the ball hits below 
the dot-dashed line in Fig. 2(a)), and the layer stays 
on the container until the acceleration becomes —g again 
(the intersection of the dot-dashed line and the trajec- 
tory of the plate). The layer leaves the container at the 
same phase angle of the oscillation at every cycle; hence 
the take-off velocity or the flight time is single- valued. 
This regime is called Period 1, n = 1^ which means the 
period of the trajectory is single-valued [Period 1) and 
the ball collides with the plate every cycle [n = 1). For 
F > 4.0, the trajectory consists of two different flight 
times, and the flight time is still a fraction of the period 
T (Fig. 2(b), Period 2, n = 1); it corresponds to the //2 
hexagonal pattern. In this regime, the magnitude of the 
container acceleration at collision is larger than g once 
every other cycle (the ball hits above the dot-dashed line 
in Fig. 2(b)). At this collision, the layer leaves the con- 
tainer immediately, and the take-off velocity and flight 
time are smaller than those of the other cycles, in which 
the ball stays on the container until the acceleration of 
the container becomes —g. 

For values of F above 4.5, the flight time of the layer 
exceeds T, and the layer hits the container once every 
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other cycle. In this regime, the layer can hit the con- 
tainer either on odd cycles or on even cycles; the tra- 
jectory becomes degenerate (the trajectory in the //2 
hexagonal pattern regime is also degenerate, and a //2 
hexagonal pattern can have a phase discontinuity defect 
which is different from a kink). The trajectory of sin- 
gle inelastic ball in the //4 square/stripe pattern regime 
(r > 5.5) is shown in Fig. 2(c) {Period 1, n = 2). When 
r is increased to 6.5, the trajectory of the layer consists 
of longer and shorter flight times, and a bifurcation to 
the //4 hexagonal pattern occurs (Fig. 2(d), Period 2, 
n = 2). At r ~ 8.0, another bifurcation to Period 1 
and n = 3 state occurs, which corresponds to a //3 flat 
pattern; however, the //3 flat pattern has not been ob- 
served before, and the layer has been known to exhibit 
spatiotemporal chaos in this regime. We found transient 
//3 and //6 patterns in this study (see Sec. VI). 

Since the layer collides with the container every other 
cycle for F > 4.5, domains tt out of phase may coexist 
in the layer. When these opposite phase domains co- 
exist, there is a phase discontinuity line defect between 
the adjacent phase domains, which we call a kink fl^ ] 
(two sides of a phase discontinuity defect in the //2 
hexagonal pattern regime are not at opposite phases, 
and this defect is not a kink). In the experiments, kinks 
are created by inhomogeneous initial conditions or ex- 
ternal perturbations such as side wall friction or tilt of 
the container; see Sec. IV. A sequence of an //4 pattern 
with kinks obtained from a 2D simulation is shown in 
Fig. 3, where one domain has fully developed pattern, 
and the other is nearly flat. An animated movie is in 
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B. Phase bubbles and randomly moving labyrinths 

As F is increased further from the / /4 hexagonal pat- 
tern regime, the layer exhibits spatiotemporal chaos, 
which is not included in the single ball model. Snap- 
shots obtained from the experiments and simulations of 
this regime are shown in Fig. 4. 

As F is increased above a critical value Fpf, in the //4 
hexagonal pattern regime, a small localized region spon- 
taneously changes its phase angle, and this region be- 
comes enclosed by a kink. We call this localized region 
a phase bubble; see Fig. 4(c) and 4(d). Phase bubbles 
suddenly (in one cycle) pop up at random locations, and 
they shrink and disappear over several tens to hundreds 
of cycles. The nucleation rate increases with F, and the 
decay time depends on the control parameters and the 
initial size. Fpf, varies with the depth of the layer and 
the material used; it is 7.5 in the experiment with bronze 
particles and N = 5, and is 7.1 in layers of lead particles 
with A'^ = 5, in 3D. In the simulations we find that the 
value of Fpf, is a few percent smaller in 2D or quasi-2D 
layers than in 3D. A phase bubble is not just a defect in 
the //4 hexagonal pattern of this system; it arises due 



to an inherent instability of the oscillated granular layer 
(see Sec. IV). Phase bubbles are spontaneously nucle- 
ated even without side wall friction or any other external 
perturbation, while kinks for 4.5 < F < Tpb are never 
created spontaneously without an external perturbation. 

As F is further increased, the nucleation rate of phase 
bubbles increases faster than the decay rate, and the layer 
eventually exhibits spatiotemporal chaos in the form of 
randomly moving labyrinths (Fig. 4(e) and 4(f)). Note 
that there are a few kinks of closed form; i.e., phase 
bubbles. In this regime, labyrinthian kinks move around 
in the container in a seemingly random fashion. Each 
phase domain of labyrinths collides with the container 
every other cycle, and an //4 subharmonic pattern is su- 
perposed on it. The layer behavior does not show any 
qualitative difference up to F ~ 14, the highest value 
investigated in the experiment. 

As a phase bubble shrinks, it becomes more circu- 
lar, as though it had a surface tension. This behav- 
ior is more clearly observed with phase bubbles in the 
//2 flat pattern regime, because there is no superposed 
pattern. Phase bubbles are not nucleated spontaneously 
for F < Tpb, but we can make them from initial condi- 
tions with inhomogeneous phase, by suddenly decreasing 
F from the randomly moving labyrinths regime to the 
//2 flat pattern regime. Such phase bubbles in the //2 
flat pattern regime are shown in Fig. 5. An animated 
movie is in 
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A kink acts on the pattern as if it were a boundary, 
and the pattern tries to rotate perpendicular to the kink. 
These "boundaries" of the layer have irregular shape; 
a well-ordered hexagonal pattern does not form in the 
phase bubble regime because phase bubbles break the 
long range order of the pattern (Fig. 4(c) and 4(d)). 
Pattern selection is also affected by the size of a phase 
bubble, or the width of a phase domain in the randomly 
moving labyrinths regime, because of the "boundary con- 
dition" of the pattern imposed by the kink. When a 
phase bubble is larger than the wavelength of the hexag- 
onal pattern, the phase bubble has an //4 hexagonal 
pattern superposed on it. Otherwise, an //4 stripe pat- 
tern is superposed because it is the only pattern that can 
fit in such small domains. 



IV. NUCLEATION OF PHASE BUBBLES 

We find in the simulation that for F slightly below Fp;,, 
the bottom of the layer exhibits an undulation, of which 
length scale is much larger than the wavelength of the 
pattern. In this Section, we discuss how the undulation 
leads to the nucleation of phase bubbles. 

In the //4 hexagonal pattern and the phase bubble 
regime, the temporal dynamics of the layer is sensitive 
to a small change in the flight time. This is the case also 
in the single ball model; in the //4 hexagonal pattern 
regime, the ball immediately leaves the container at ev- 



3 



cry other collision, so that the instantaneous velocity of 
the container at this collision becomes the take-off veloc- 
ity of the ball. Thus, if the velocity is changed slightly 
and the flight time is increased, the take-off velocity at 
the next flight becomes smaller; this velocity may be too 
small for the ball to fly over the container during the 
following two cycles, so that the ball collides with the 
container in the very next cycle (Fig. 6). In a real layer, 
the flight time has some fluctuation due to the undula- 
tion; the mechanism in Fig. 6 will create a phase bubble 
when the undulation is big enough (Fig. 7). 

The power spectral densities (PSD's) of the top and 
the bottom of a layer are obtained from the simulation. 
Spectra for the values of T slightly below Fpf, show that 
another peak with small wave number appears as F ap- 
proaches Fpb (Fig. 8). For smaller F (F < 6.5 for this 
case), the shape of the bottom of the layer is slaved to the 
top of the layer and has the same wavelength as the top; 
PSD's have peaks only at the wavelenght of the pattern 
and its subharmonics. 

The nuclcation rate of phase bubbles monotonically in- 
creases with F; we hypothesize that the amplitude or the 
growth rate of the undulation increases with F. Above 
some value of F (> Fpb), the nucleation rate exceeds the 
decay rate, and phase bubbles or kinks accumulate in the 
layer until each phase domain reaches the shortest length 
scale that cannot have another kink in it; i. e., this sets 
the width of the labyrinths. As a result, the layer is filled 
with a superposition of phase bubbles which connect to 
one another, forming non-closed kinks. We propose that 
the randomly moving labyrinth pattern results from the 
above mechanism, i.e., randomly moving labyrinths may 
be understood as a "saturated" state of phase bubbles. 

There are two more mechanisms responsible for the 

creation of phase bubbles or kinks: (1) side wall friction 
and (2) tilt of the layer. In an infinitely extended layer 
without side walls, kinks are formed only by the undu- 
lation of the layer, or nucleation of phase bubbles, if the 
layer is perfectly level. In the simulation, side walls are 
easily eliminated by imposing horizontal periodic bound- 
ary condition; however, the layer size is always finite in 
the experiment, and the side wall effect cannot be com- 
pletely eliminated. In the experiment, the friction due to 
the side walls disturbs the oscillatory motion of the layer, 
which often creates kinks. This effect becomes more and 
more important as F is increased. There arc some phase 
bubbles in contact with the side wall in Fig. 4(c) and 
4(d), which are kinks created by the side wall friction. 
In the randomly moving labyrinths regime, the side wall 
effect is one of the major reasons for the decay of //3 
or //6 transient patterns (see Sec. VI). In addition to 
the above effect, tilt of the container also creates kinks. 
This is the main mechanism of the creation of kinks for 
4.5 < F < Tpb in the experiment. 



V. DYNAMICS OF A PHASE BUBBLE 

Once formed, a phase bubble shrinks as if it had a 
surface tension, and then disappears (Fig. 5). We now 
discuss why this happens. 

When the layer collides with the plate, a density wave 
forms and propagates across the boundary of a phase 
bubble (a kink). This density wave initiates a coUisional 
momentum transfer in the horizontal as well as in the 
vertical direction, across a kink (Fig. 9). A sequence of 
side views of a kink in a 2D layer is in Fig. 10, which 
shows the relation for the density wave to the momentum 
transfer. At i = in Fig. 10, the left half of the layer 
collides with the plate (horizontal bar) and is pushed up. 
It becomes compact and nearly static with respect to the 
plate, acting as if it were a part of the container. At this 
moment, the rest of the layer is still falling and is dilated. 
The two parts of the layer interact at the kink, and a den- 
sity wave forms at the interface due to the large density 
gradient. The density wave propagates toward the di- 
lated part until the compact part loses contact with the 
plate. The density wave drives a coUisional momentum 
transfer and pushes the interface toward the dilated part 
(the kink is shifted rightward at i = 1.2T compared to 
its earlier location at i = 0) . In the next cycle, the same 
process occurs in the opposite direction, and the kink is 
pushed back to its original position (not shown in Fig. 
10), provided that the phase difference is exactly tt (then 
both parts of a 2D layer are symmetric). As a result, 
the momentum flux across a kink in a 2D or quasi-2D 
layer oscillates symmetrically with a period of 2T, with 
no net translational motion; thus a phase bubble in a 2D 
or quasi-2D layer never shrinks but only oscillates sym- 
metrically (Fig. 11(a)). 

The mechanism of the density wave and the momen- 
tum transfer across a kink in a 3D layer is the same as 
in a 2D/quasi-2D layer; however, in a 3D layer, there 
is a geometric effect which is absent in a 2D/quasi-2D 
layer: the two sides of a kink are not symmetric unless 
the kink is straight. As the particles at the front of the 
density wave are nearly static with respect to the plate, 
the momentum fiux is roughly proportional to the num- 
ber of particles at the front of the density wave and the 
velocity of the container. For a kink with nonzero local 
curvature, the number of particles at the front of the den- 
sity wave for the two successive cycles are different, even 
if the two sides of the kink are exactly tt out of phase. 
Thus the momentum transfer in one direction is always 
larger than that of the other. This leads to an asym- 
metric momentum transfer, or an asymmetric oscillatory 
motion of a kink. Consequently, the total net momen- 
tum flux across the boundary of a phase bubble over a 
multiple of 2T is always inward, because its boundary is 
a kink of a closed shape. As a result, a phase bubble 
shrinks and disappears. An asymmetric oscillation of a 
boundary of a phase bubble is shown in Fig. 11(b). The 
difference in the numbers of particles at the interface of 
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the density wave for the two consecutive cycles is approx- 
imately proportional to the local curvature of the kink. 
Thus, the speed of a kink is roughly proportional to the 
local curvature, leading to an effective surface tension. 

There are other mechanisms driving the motion of a 
kink in the experiment, which are minor compared to 
the previously discussed curvature-dependent geometri- 
cal effect: (1) non-iT phase difference effect and (2) the 
finite mass ratio effect. Firstly, we have assumed in the 
above discussion that both sides of a kink are exactly tt 
out of phase; if this were not the case, the impact from 
the container for two successive cycles would be different, 
and the oscillation of a kink would be asymmetric, even 
if the local curvature were zero. This is why a secondary 
forcing can control the motion of a straight kink, as was 
experimentally observed by Aranson et al. . Secondly, 
if the difference in the mass of each phase domain of the 
layer is not negligible compared to mass of the container, 
the impact by the container on both domains of the layer 
are different, and a net translational motion of a kink is 
induced, even if it has zero local curvature and the phase 
difference is tt. We call this the finite mass ratio effect. 
In the experiment, a kink in the //2 flat pattern regime 
eventually travels to the center of the container, until the 
impacts from the container on both phase domains are 
balanced; it was observed in Refs. jl^Jl^. This effect is 
absent in the simulation, because the mass of the con- 
tainer is assumed to be infinitely large compared to the 
mass of the layer. 



VI. TRANSIENT //3 AND //6 SUBHARMONIC 
PATTERNS 



Until now, no standing wave pattern has been observed 
above the //4 hexagonal pattern regime in the phase di- 
agram in Fig. f . In this Section we present the first ob- 
servation of the //3 flat and f/6 square/stripe patterns. 
These patterns were found first in the simulations, which 
motivated the experiments. We start the discussion with 
the single ball model. 

The single ball model predicts an infinite cascade of 
bifurcations with increasing natural number n, for in- 
creasing r (there are small chaotic windows along F, but 
this is not important in this discussion). For the granu- 
lar layers, these bifurcations correspond to the following 
patterns: 



f/n flat ^ 

1 

//(2n) squares/stripes 

1 

f/{2n) hexagons 



n ^ n + 1 



where the initial n is I. This model predicts a bifurca- 



tion from an //4 to an //3 state {Period 1, n = 3) dX 
around F = 8.0, however, in the experiments, the cascade 
of bifurcation stops at n = 2. 

We discussed how the undulation of the layer leads to 
the nucleation of phase bubbles in Sec. IV. If the un- 
dulation of the layer could be avoided, we expect that 
//3 flat and //6 square/stripe patterns would exist. To 
test this conjecture, we performed a series of numerical 
simulations of a small layer (IOct x lOcr and 20(7 x 20(7 
with A'^ = 6) in the randomly moving labyrinths regime 
for (F, /*) = (8 - II, 0.8 - 1.0), using periodic horizontal 
boundary conditions. In these simulations, we observed 
that the layer collides with the plate every third cycle, 
and that the dynamics is stable up to order of 1000 cycles; 
we conclude that if a larger layer followed this behavior, 
//3 or //6 subharmonic patterns would exist. 

In larger layers, the undulation of the layer cannot be 
completely avoided, but can be suppressed at least dur- 
ing a short time by preparing a flat and compact layer as 
an initial condition; when F is quickly increased from the 
flat pattern regime to where //3 or //6 subharmonic pat- 
terns are predicted by the single ball model (8 ;$ F 11), 
a //3 or //6 pattern is found in the simulation. Later 
these transient patterns were found in the laboratory ex- 
periments as well, as Fig. 12 and 13 illustrate. These 
patterns emerge as the primary instability of the layer in 
this regime, but they are invaded and overtaken by either 
kinks formed due to the side wall friction (Fig. 12) or by 
the undulation of the layer (Fig. 13). As a result, the do- 
main of these transient patterns gradually decreases and 
is eventually taken over by randomly moving labyrinths. 
These transient patterns last up to several hundreds of 
cycles, depending on the initial condition, control param- 
eters, and the system size. 



VII. DISCUSSION 

We have shown that phase bubbles play a critical role 
in the transition to spatiotemporal chaos in the patterns 
formed by vertically oscillated granular layers. Phase 
bubbles spontaneously form first in the f/4 hexagonal 
pattern regime as the acceleration F is increased. At 
larger F, the rate of nucleation of bubbles grows faster 
than the decay rate, leading ultimately to a spatiotem- 
porally chaotic pattern of randomly moving labyrinths. 
In a qualitatively similar way, the formation of defects 
has been found to lead to spiral defect chaos Q and to 
chaos in a model with invasive defects ||l^ . 

We have investigated the mechanism of the nucleation 
and the dynamics of phase bubbles and randomly moving 
labyrinths, using inelastic hard sphere molecular dynam- 
ics simulations and experiments. We have found that a 
vertically oscillated granular layer has a large scale un- 
dulation, even without interstitial gas (Sec. IV). The 
undulation constitutes an inherent feature of oscillated 
granular layers, like the standing wave pattern forma- 
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tion, and cannot be avoided; above some critical value 
of r in the //4 hexagonal pattern regime, the undula- 
tion of the layer gives rise to the nucleation of a phase 
bubble. The spontaneous nucleation of kinks in verti- 
cally oscillated 2D or quasi-2D granular layers for large 
r has been called "subharmonic instability" ||l9| or "arch- 
ing" pO| , pl| , but the mechanism has not been elucidated. 
These two phenomena are all essentially phase bubbles 
in 2D or quasi-2D oscillated granular layers. 

We also found that a kink of nonzero curvature has 
a net translational motion due to the asymmetric colli- 
sional momentum transfer across a kink (Sec. V). The 
local speed of a kink is roughly proportional to the local 
curvature; a phase bubble shrinks as if it had a surface 
tension and then disappears, because it is a kink of a 
closed form. We showed that shrinking of a phase bub- 
ble and translational motion of a kink are essentially the 
same phenomena. 

Finally, based on the understanding of kinks, we pre- 
dicted transient //3 and //6 subharmonic patterns and 
observed them for the first time (Sec. VI). 
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FIG. 1. Phase diagram obtained from an experiment with 
bronze particles of diameter cr = 165 /um and nondimensional 

depth A'^ = 5, in a circular container with diameter L = 770a, 
showing particularly the details for F > 7.0. Solid lines de- 
note the transitions for increasing F, while the dotted lines 
denote decreasing F. 
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a) r = 3.0, f/2 square/stripe regime 




1 2 3 

(b) r = 4.0, f/2 hexagon regime 




1 2 3 

(c) r = 6.0, f/4 square/stripe regime 
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(d) r = 7.0, f/4 liexagon regime 




4 6 8 

time (t/T) 

FIG. 2. Temporal trajectory of a completely inelastic ball 
(e = 0) on an oscillating plate, which models the trajectory of 
the center of mass of the layer. The solid sinusoidal curve is 
the trajectory of the plate. The ball leaves the plate when the 
acceleration of the plate becomes — g; i. e., where the hor- 
izontal dot-dashed line intersects with the trajectory of the 
ball. If the ball collides with the plate above the dot-dashed 
line (in (b) and (d)), it leaves the plate immediately. 
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lOOo 200(5 



FIG. 3. A side view of a 2D layer of an //4 pattern with 
two kinks, obtained in the simulation for F = 6.5, /* = 0.8, 
and N = 8. The middle part and the rest of the layer oscillate 
TT out of phase, and phase discontinuities between these two 
phase domains arc called kinks (indicated by arrows). The 
container (horizontal bar at the bottom) is at its minimum 
height at the phase angle at which these figures are taken. 
Horizontal boundary is periodic. 
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EXPERIMENTS SIMULATIONS 

f/4 hexagonal pattern 







Phase bubbles in f/4 hexagonal pattern 




Randomly moving labyrinths 




FIG. 4. Top views obtained from the experiments, and 
the mean height as a function of the position obtained from 
3D simulations. In (a) and (c), there is a phase kink near 
the middle, while there is not in (b) or (d). Kinks are eas- 
ily created in the experiment due to the side wall friction or 
tilt of the layer, but they are easily eliminated in the simula- 
tion; see Sec. IV. One phase bubble is indicated by a white 
arrow in (c); in (d) the phase bubbles are the white areas, be- 
cause they reach their maximum height at this moment and 
patterns superposed on them are nearly flat. The gray scale 
in the experimental images indicates the intensity of the re- 
flected light, which is a measure of the gradient of the surface. 
The gray scale in the simulation images is {z){x, y), increasing 
from black to white, where (F){x, y) is an averaged value of F 
over the particles located at {x, y); {z){x, y) is called the mean 
height field. A circular container of diameter L = 847a is used 
in the experiments, and a circular container of L = 294(T is 
simulated. The parameters (r,/*,iV) are (a) (7.1,1.0,10), (b) 
(7.0,0.85,8), (c) (7.3,1.0,10), (d) (7.2,0.85,8), (e) (8.0,1.2,15), 
and (f) (8.9,0.85,8). 
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t = S"*^^" 




FIG. 5. A sequence of phase bubbles in //2 flat pattern, 
obtained in the experiment. The bubbles shrink as if they 
had a surface tension. To create these phase bubbles, F was 

suddenly decreased from 9.0 to the //2 flat pattern regime, 
F = 4.5 (/* = 0.6). Experimental setup is the same as in Fig. 
1. 400cr X 400cr portion of the layer is shown. 
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time (tjT) 

FIG. 6. Sensitivity of the trajectory to small perturbar- 
tions is illustrated by perturbed (dashed line) and unper- 
turbed (solid line) trajectories of a completely inelastic ball 
in the //4 hexagonal pattern regime (F = 7.0). For the per- 
turbed trajectory, we increased the initial take-off velocity by 
3% to delay the collision. At the next collision (indicated by 
an arrow), the perturbed ball collided with the container be- 
cause its take-off velocity was too small, and its trajectory 
became out of phase with that of unperturbed one. 
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FIG. 7. (a) A sequence of side views of a quasi-2D layer 
showing the nucleation of a phase bubble around B, obtained 
in a simulation of a quasi-2D layer of size 200(7 x lOcr, for 
r = 7.2, /* = 0.9, and N = 6. (b) The trajectories of the 
center of mass at two horizontal locations, A and B in (a), 
as a function of time. The sinusoidal function shown in (b) is 
the trajectory of the plate of the container. At t = in (a), 
the bottom of the layer has a short length scale deformation 
slaved to the pattern on the surface of the layer. This undula- 
tion at the bottom grows, as shown in the next two successive 
frames in (a) {t = 4T and GT in (b)). At t = 6T, the colhsion 
of B is delayed due to the undulation; at this collision, a por- 
tion of the layer around B collides with the container when 
the plate nearly reaches its maximum height, and the take-off 
velocity becomes too small to fly the next two cycles over the 
container. As a result, B collides with the container at the 
very next cycle, at f ~ 7T, and it becomes nearly tt out of 
phase with the rest of the layer. This mechanism is the same 
as that in Fig. 6. 
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FIG. 8. Power spectral densities (PSD's) of the top and 
the bottom of a quasi- 2D layer for F = 6.6, /* = 0.9, and 
N = 6, obtained in simulations of a quasi-2D layer of size 
lOOOo" X 10a having more than 50 wavelengths; there is an 
additional peak (indicated by an arrow) at a small wavenum- 
ber as well as the primary peaks at the wavenumber of the 
pattern (at ko — l-najX, where A is the wavelength of the 
pattern) and its subharmonic. Fpb = 6.7 for this simulation. 
PSD's were measured when the amplitude of the pattern weis 
fully developed, just before the layer collided with the con- 
tainer. 
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FIG. 9. (a) Horizontal momentum field {\/ pi -\- Py){x,y), 
which has its maximum values (indicated white) along the 
boundary of phase bubbles; there is a significant momentum 
flux perpendicular to the boundary (the kink), (b) Vertical 
momentum field {pz){x,y), black being the maximum down- 
ward, and white being the maximum upward. Fluctuation in 
vertical velocities (blackish areas) is due to the large scale un- 
dulation of the layer, discussed in Sec. IV. These are obtained 
from the same simulation data as in Fig. 4(d). 
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FIG. 10. A sequence of an //2 flat pattern with a kink, ob- 
tained from the simulation of a 2D layer for T = 5.2, /* — 0.6, 
and N = 9. Each circle corresponds to a particle. The gray 
scale of the particles indicates the magnitude of the horizontal 
momentum \px\, increasing from black to white. As the layer 
is pushed up by the plate, a density wave forms, which ini- 
tiates a momentum transfer propagating rightward (see text 
for details). 
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X (a) X (a) 

FIG. 11. Space-time plots of horizontal momentum 
(Px){x,t) of a phase bubble (a) in a quasi-2D f/2 flat pat- 
tern, and (b) in a 3D f/2 flat pattern (cross section along the 
diameter of a phase bubble). The values increase from black 
to white, black being the maximum leftward and white being 
the maximum rightward. A phase bubble in a quasi- 2D layer 
oscillates symmetrically, while in a 3D layer oscillates asym- 
metrically, shrink and disappear (indicated by an arrow). A 
quasi- 2D layer of size IOOct x 10a is used for (a), and the 
phase bubble in (b) was created by decreasing F from 8.5 to 
5.2. /* = 0.6 and AT = 6 for both. 
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t = 72T 



t = 1 03T 




FIG. 12. A sequence of a transient //3 flat pattern ob- 
tained from the experiment. The friction due to the side wall 
creates a kink, which propagates to the center of the container 
in the radial direction (the same mechanism of the shrinking 
of a phase bubble) and destroys the //3 flat pattern. T was 
suddenly increased from below the onset, about 2.0, to 7.8; 
/* = 0.94 and A'^ = 6. The experiment was done with the 
same particles as in Fig. 1, and L = 847cr. 



t = 1 = 1 2T 




FIG. 13. A sequence of a transient //6 pattern, obtained 

from the experiment. The undulation of the layer creates 
kinks in the middle of tiie layer, which destroys the //6 pat- 
tern, r was suddenly increased from 2.0 (below pattern on- 
set), to 10 (/* = 1.2 and N = 6). The experiment was done 
with the same particles, N, and L as in Fig. 13, and a part 
of the layer of size 250ct x 250cr is shown. 
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